setwd("your/local/directory/")

library("foreign")

ees <- as.data.frame(read.dta("EES 2009 - stacked.dta"))
attach(ees)

stderr <- function(x) sqrt(var(x)/length(x))
####################################
postscript("fig4.eps",width=10,height=3.75)
par(oma=c(3.1,3.1,0,0))
par(mar=c(1, 1, 1, 1))
par(cex=.4)
par(mfrow=c(2,5))
data <- list()
avg <- list()
avg.se <- list()
avg.hi <- list()
avg.lo <- list()
smple <- list()
ylims <- c(0,8)
for(i in c(-1:-5,1:5)+6) {
  data[[i]] <- as.data.frame(subset(ees, lrcen == (i-6) & complete.cases(ptv,stack,lrcen) & eu15 ==1))
  smple[[i]] <- sample(1:NROW(data[[i]]), 1000, replace = FALSE)
  avg[[i]] <- tapply(data[[i]]$ptv,data[[i]]$lrprtycen,mean)
  avg.se[[i]] <- tapply(data[[i]]$ptv,data[[i]]$lrprtycen,stderr)
  avg.hi[[i]] <- avg[[i]] + avg.se[[i]]*1.96
  avg.lo[[i]] <- avg[[i]] - avg.se[[i]]*1.96
  avg[[i]] <- c(avg[[i]][1],avg[[i]][1:11])
  avg.hi[[i]] <- c(avg.hi[[i]][1],avg.hi[[i]][1:11])
  avg.lo[[i]] <- c(avg.lo[[i]][1],avg.lo[[i]][1:11])
  plot(data[[i]]$lrprtycen[smple[[i]]],data[[i]]$ptv[smple[[i]]], type="n",ylim=ylims, xlab="Left-Right Party-Placement", ylab="Party Preference Strength", axes=F)
  polygon(c(c((i-6)-.5,(i-6)-.5),c((i-6)+.5,(i-6)+.5)),c(c(ylims[2],ylims[1]), c(ylims[1],ylims[2])),
        col = "grey90", border = NA)
  abline(v=-.5, lty=2, lwd=1, col="darkgrey")
  abline(v=.5, lty=2, lwd=1, col="darkgrey")
  lines(-5.5:5.5,avg[[i]], lwd=1,type="S")
  lines(-5.5:5.5,avg.lo[[i]], lwd=.7,type="S",lty=2)
  lines(-5.5:5.5,avg.hi[[i]], lwd=.7,type="S",lty=2)
  if (i == 5 | i == 7) {axis(2)}
  if (i > 6) {axis(1, at=c(-5,-3,-1,1,3,5), labels=c(-5,-3,-1,1,3,5))}
}
mtext("Left-Right Party-Placement", side=1, outer=TRUE, adj=.5,cex=.8,line=1.5)
mtext("Party Preference Strength", side=2, outer=TRUE, adj=.5,cex=.8,line=1.5)
dev.off()
####################################
postscript("fig3.eps",width=10,height=3.75)
par(oma=c(3.1,3.1,0,0))
par(mar=c(1, 1, 1, 1))
par(cex=.4)
par(mfrow=c(2,5))
data <- list()
avg <- list()
avg.se <- list()
avg.hi <- list()
avg.lo <- list()
smple <- list()
ylims <- c(0,8)
for(i in c(-1:-5,1:5)+6) {
  data[[i]] <- as.data.frame(subset(ees, lrprtycen == (i-6) & complete.cases(ptv,stack,lrcen,lrprtycen) & eu15 ==1))
  smple[[i]] <- sample(1:NROW(data[[i]]), 1000, replace = FALSE)
  avg[[i]] <- tapply(data[[i]]$ptv,data[[i]]$lrcen,mean)
  avg.se[[i]] <- tapply(data[[i]]$ptv,data[[i]]$lrcen,stderr)
  avg.hi[[i]] <- avg[[i]] + avg.se[[i]]*1.96
  avg.lo[[i]] <- avg[[i]] - avg.se[[i]]*1.96
  avg[[i]] <- c(avg[[i]][1],avg[[i]][1:11])
  avg.hi[[i]] <- c(avg.hi[[i]][1],avg.hi[[i]][1:11])
  avg.lo[[i]] <- c(avg.lo[[i]][1],avg.lo[[i]][1:11])
  plot(data[[i]]$lrcen[smple[[i]]],data[[i]]$ptv[smple[[i]]], type="n",ylim=ylims, xlab="Left-Right Self-Placement", ylab="Party Preference Strength", axes=F) 
  polygon(c(c((i-6)-.5,(i-6)-.5),c((i-6)+.5,(i-6)+.5)),c(c(ylims[2],ylims[1]), c(ylims[1],ylims[2])),
          col = "grey90", border = NA) 
  abline(v=-.5, lty=2, lwd=1, col="darkgrey")
  abline(v=.5, lty=2, lwd=1, col="darkgrey")
  lines(-5.5:5.5,avg[[i]], lwd=1,type="S")
  lines(-5.5:5.5,avg.lo[[i]], lwd=.7,type="S",lty=2)
  lines(-5.5:5.5,avg.hi[[i]], lwd=.7,type="S",lty=2)
  if (i == 5 | i == 7) {axis(2)}
  if (i > 6) {axis(1, at=c(-5,-3,-1,1,3,5), labels=c(-5,-3,-1,1,3,5))}
}
mtext("Left-Right Self-Placement", side=1, outer=TRUE, adj=.5,cex=.8,line=1.5)
mtext("Party Preference Strength", side=2, outer=TRUE, adj=.5,cex=.8,line=1.5)
dev.off()
